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Q ' Abstract. Several surfactant molecules self-assemble in solution to form long, flexi- 

. ble wormlike micelles which get entangled with each other, leading to viscoelastic gel 

phases. We discuss our recent work on the rheology of such a gel formed in the dilute 
aqueous solutions of a surfactant CTAT. In the linear rheology regime, the storage 
modulus G'{lu) and loss modulus G"{uj) have been measured over a wide frequency 
'■^ , range. In the nonlinear regime, the shear stress a shows a plateau as a function of 

Ch ' the shear rate 7 above a certain cutoff shear rate 7c. Under controlled shear rate 

conditions in the plateau regime, the shear stress and the first normal stress difference 
show oscillatory time-dependence. The analysis of the measured time series of shear 
stress and normal stress has been done using several methods incorporating state space 
reconstruction by embedding of time delay vectors. The analysis shows the existence 
of a finite correlation dimension and a positive Lyapunov exponent, unambiguously 
' implying that the dynamics of the observed mechanical instability can be described by 

, that of a dynamical system with a strange attractor of dimension varying from 2.4 to 

00 ■ 2.9. 
O 
O 
O 



I INTRODUCTION 

The amphiphilic nature of surfactant molecules results in their reversible self- 
. assembling to form larger structures like spherical and cylindrical micelles, bi- 

Q \ layers and vesicles. These supramolecular assemblies self-organise to form either 

J-^ ' long-range ordered liquid crystalline phases or short-range ordered isotropic liquid 

phases. Some examples of ordered phases are the nematic ordering of micellar discs 
^ ! or rods, the smectic A phase formed from the one- dimensional packing of infinite 

^ \ lamellae, cubic phases formed by finite-size spherical micelles and two-dimensional 

hexagonal crystals formed by infinite cylindrical micelles. In recent years, bicontin- 
uous cubic phases have been observed in which a continuous bilayer or a monolayer 
forms a periodic, saddle-shaped, minimum energy surface. In isotropic phases, mi- 
croemulsions and sponge phases can form. Another example is the formation of 
viscoelastic gels from the entanglement of very long and semiflexible cylindrical 
micelles with polymerlike behaviour. Interestingly, these micelles break and rejoin, 
giving rise to the name 'wormlike micelles'. The length distribution depends on 



the surfactant and the salt concentrations, temperature and the energy of scission 
of the micelles. The objective of this paper is to discuss our recent results on the 
nonlinear flow behaviour of such a system formed by the surfactant CTAT (cetyl 
trimethylammonium tosilate) in water. At very low concentrations, (c < 0.04wt.%) 
and above the Kraft temperature of 23° C, CTAT molecules in water self- assemble 
to form spherical micelles exhibiting Newtonian flow. Above 0.9wt.%, CTAT forms 
giant wormlike micelles. The hexagonal phase is seen at a CTAT concentration of 
above 27wt.% [1,2]. 

In polymeric melts, stress relaxation occurs primarily by reptation dynamics 
(curvilinear diffusion of the polymeric chain along the contour of an imaginary 
tube enclosing the chain), as suggested by de Gennes [3], which defines a relaxation 
time Trep. In systems of giant wormlike micelles, there is an additional contribution 
to the stress relaxation from the reversible scission (reversible breakage and recom- 
bination) of micelles, which occurs on a time scale r^. If t;, << r^ep, then stress 
relaxation in wormlike micelles occurs mainly by reversible scission and is described 
well by the Maxwell model which gives the response function G(t)=Goea;p(— ;^), 

where Go is the high frequency plateau modulus and tm — {TbTrepY^'^ is the relax- 
ation time. The physical picture of the single relaxation time is a kind of motional 

narrowing [4,5]. 

The reptation reaction model of wormlike micelles, which gives the Maxwell 
model for the linear viscoelasticity, yields unusual nonlinear rheology observed in 
some viscoelastic gels (for example cetyl pyridinium chloride/ sodium salicylate) 
[6,7]. The flow curve is characterised by a plateau in the shear stress versus shear 
rate curve [6], while the normal stress is found to increase linearly with shear 
rate [6,8]. These features in the flow curve of giant wormlike micellar systems are 
associated with a flow instability of the shear banding type, where the sheared 
solution breaks up into bands of high and low viscosities supporting low and high 
shear rates, respectively [8] . Flow birefringence [9] and nuclear magnetic resonance 
velocity imaging [10] have revealed the existence of banded flow in the shear stress 
plateau regime. An alternate hypothesis explains the plateau in the shear stress 
as due to an isotropic-nematic phase transition that occurs due to the application 
of high shear rates to the sample [11]. Experiments on the transient rheology of 
CPyCl-NaCl after imposition of step shear rates, show overshoots, oscillations and 
long time sigmoidal kinetics in the shear stress a{t) [11]. In our system of interest, 
CTAT 1.9wt.% in water at 25°C, we do not observe any overshoots in the stress 
relaxation on application of shear rates lying in the plateau region of the flow 
curve. However, we see that at high shear rates, the stress instead of decaying to 
a steady state, oscillates in time for a long time. We have analysed the time series 
of the viscoelastic and normal stresses obtained from these experiments and have 
demonstrated the existence of deterministic, chaotic dynamics in sheared solutions 
of CTAT. Oscillations in the stress relaxation of CTAT have also been observed by 
Soltero et al at a concentration of 5wt.% [12]. 

Chaotic time series data have often been observed in experiments on fluid dy- 



namics [13]. The Couette-Taylor flow in fluids and the Rayleigh-Benard convection 
in gases are two examples. Ananthakrishna et al [14] have shown the existence 
of chaotic dynamics in the jerky flow (Portevin-Le Chateher effect) of some metal 
alloys undergoing plastic deformation. These systems also exhibit nonmonotonic 
flow curves. Deterministic chaos represents the apparently irregular behaviour of 
dynamical systems that arises from strictly deterministic laws in the absence of any 
external stochasticity. Chaotic dynamics in physical systems is characterised by an 
exponentially sensitive dependence on initial conditions, as a result of which long- 
term predictability of the dynamics of these systems is impossible. Chaotic time 
series may be characterised by certain invariants, metric and dynamical, such as 
the various fractal dimensions, the largest Lyapunov exponent and the Kolmogorov 
entropy. A positive value of the largest Lyapunov exponent is a direct consequence 
of the sensitivity of the trajectories in phase space to small changes in the initial 
conditions. The fractal dimensions provide a measurement of the topology of the 
attractor on which all the trajectories asymptotically converge. Apart from describ- 
ing the complexity of the attractor, they also provide a measure of the regularity 
of occurrence of the phase points on its surface [15]. 

In this paper, we will show the existence of a positive Lyapunov exponent and a 
flnite correlation dimension in the observed time series of the stress relaxation of 
CTAT on application of step shear rates lying in the plateau region of the flow curve. 
Further, the magnitudes of the invariants are found to increase monotonically with 
shear rate, the control parameter in our experiments. 



II EXPERIMENTAL DETAILS 

CTAT (purchased from Sigma Chemicals, India) samples were prepared by dis- 
solving 1.9wt.% of powdered CTAT in distilled and deionised water. The samples 
were kept in an incubator at a temperature of 60°C for a week, and frequently 
shaken to facilitate homogenisation. The samples were then kept at the experi- 
mental temperature (25° C) for two days. The rheological properties of the CTAT 
thus prepared and equilibrated were measured using a Rheolyst AR-IOOON ( T. 
A. Instruments, U. K.) stress-controlled rheometer with temperature control and 
software for shear rate control. The rheometer was also equipped with four strain 
gauge transducers at a distance of R from the centre of the plate for the measure- 
ment of the normal force F^. is related to the flrst normal stress difference Z by 
Z = 1^ [16]. All the measurements reported in this paper were performed using a 
conc-and-plate geometry [16] of diameter 4cm and cone angle 1°59". In each flow 
experiment, the rheometer can collect upto a maximum of 1500 data points, with 
a time gap of Is or more between acquisition of successive data points. Sample his- 
tory effects were found to be important in CTAT and hence all experiments have 
been done on fresh samples from the same batch. Most of the experiments were 
performed at 25°C and a few at 35°C. 
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FIGURE 1. Frequency response of CTAT 1.9wt.% at (a) 25°C and (b) 35°C, at an oscillatory 
stress of 0.6Pa. G'{lu) is indicated by solid circles and G"{uj) by hollow squares. The solid and 
dashed lines show the fits to the Maxwell model for G'(a;) and G"(a;), respectively. 



Ill RESULTS 



Fig. 1 shows the frequency response measurements at (a) 25°C and (b) 35°C. 
The hnear regime was first ascertained by imposing stresses between 0.005 Pa to 2 
Pa, oscillating at a frequency of 0.1 Hz, and finding the window of the applied stress 
values over which the elastic modulus G'(a') and the viscous modulus G"(ci;) are 
constant. From this measurement, the oscillatory stress for the frequency response 
measurements at both the temperatures has been chosen to be 0.6 Pa. In Figs. 
1(a) and 1(b), the solid circles and the hollow squares represent G'(a;) and G"(a;), 
respectively. The dominant relaxation time, estimated from the angular frequency 
at which G'(a;) and G"{u) are equal, is found to be 17s (Is) at 25°C (35°C). The 
solid {G'{u!)) and dashed lines {G"{uj)) are the corresponding fits to the Maxwell 

model, given by G'(uj) = and G"(uj) = Sf^. As mentioned before, 

for viscoelastic gels of wormlike micelles formed by the CTAC-NaSal,CPyCl-NaSal 
and CPySal-NaSal systems [4,7], the Maxwell model works very well. However, for 
CTAT at 25°C, the fit of the data to the Maxwell model is rather poor. Interestingly, 
the fit is reasonably good for experiments done at 35°C. This may be due to the 
fact that the motional narrowing condition, << r^ep, is not satisfied at low 
temperatures. 

Fig 2 shows the viscoelastic stress and normal stress measurements when the 
stress was increased from 10 mPa to 4.7 Pa after allowing a time tm between 

acquisition of consecutive data points, and measuring the corresponding shear rate. 
These measurements give the metastable branches of the flow curve at 25°C (Fig. 
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FIGURE 2. The metastable flow curves of CTAT 1.9wt.%, showing the shear stress and normal 

stresses, respectively, as a function of shear rate at 25°C ((a) and (b)) and 35°C ((c) and (d)). 

2(a)) and 35°C (Fig. 2(c)). Figs. 2(b) and 2(d) show the measured normal stresses 
as a function of increasing shear rate at both the temperatures. At 25°C, the 
viscoelastic stress (2(a)) shows a crossover to a plateau, while the normal stress 
increases quadratically (2(b)) with increase in shear rate. The presence of a plateau 
in the viscoelastic stress at high shear rates is a signature of a mechanical instability, 
possibly of the shear banding type [7,8]. We have ruled out the possibility of an 
isotropic- ncmatic phase transition on account of the very low volume fractions of 
CTAT used in our experiments. At 35°C, the viscoelastic stress does not show a 
plateau region and the normal stress is roughly zero, which implies the absence of 
instability in sheared solutions of wormlike micelles at higher temperatures [17]. 

We have studied in great detail the relaxation of the viscoelastic stress on ap- 
plication of step shear rates lying in the plateau region of the flow curve. Fig. 3 
(a)-(d) show the stress relaxation measured for 1500 seconds, after subjecting the 
sample to step shear rates of 30, 60, 90 and 150s^^ (Figs. 3 (a)-(d)) at 25°C. The 
stress decays initially for t~100 seconds (not shown in Fig. 3), and then oscillates 
in time. This time-dependent behaviour of the stress becomes more pronounced 
at higher shear rates. Fig. 3(c) show the stress relaxation at 35°C on applica- 
tion of a step shear rate 7 = 150s~^ The oscillations seen at 25°C at 7=150s~^ 
are found to disappear completely at 35°C, leading us to ascertain that the time- 
dependent oscillations observed in the stress at 25°C are a direct consequence of the 
mechanical instability in the system. The first normal stress difference Z, which is 
measured simultaneously with the viscoelastic stress a at 7=150s~^ also shows os- 
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FIGURE 3. Shear stress relaxation measured in CTAT 1.9wt.% at 25°C on application of a step 
shear rate of (a) 30s^^, (b) 60s^^, (c) 90s~^ and (d) 150s~^. (e) shows the relaxation of shear 
stress at 35°C and 7 = 150s~^.(f) shows the relaxation of normal stress at 25°C on application 
of a shear rate of 150s~^. 



dilations, as shown in Fig. 3(f). In what follows, we describe how we have analysed 
the time-dependent data in a and Z to establish the presence of chaotic dynam- 
ics in the sheared CTAT solution. The time series analysis has been performed 
by using the method of state space construction by embedding time delay vectors 
[18]. Suppose {uj = a{jA t), j=l,1500} denotes the time series in stress a on 
application of a step shear rate. At — Is is the time elapsed between collection of 
data points. In this method, we construst m-dimensional L-delay vectors given by 
{Xi = {(Ti, (Ti+L, (Ti+2L, o-i+(m-i)L)}- {Xi:i= 1 ,2,. . . ,N- (m- 1 )L} defines a trajectory 
in m-dimensional space due to the presence of a dynamics F(Xj) : Xi Xi^i. 
In doing so, it is crucial to choose the optimal embedding dimension rrio and the 
optimal delay time Lq properly [19]. We have done this by analysing the local ex- 
ponential divergence plots and the plots of 7 = < ln{\\Xi+k — ^j+k\\/\ > 
versus L for the trajectories obtained by embedding the time series in m dimensions. 
Here < .. > denotes averaging over all pairs. Embedding the trajectories in 
higher dimensions leads to a reduction in the number of false neighbours as at m 
< nio, we are not looking at the real dynamics, but at its projection. At m> nio, 
the number of false neighbours is zero. The left panel of Fig 4 shows the local 
exponential divergence plots for 7=150s^^, where the abscissa is ln(||Xj — Xj\\) 
and the ordinate is ln(||Xj+fc — — Xj\\) for (L,m), where L= 1 and 5 

and m=3 to 6. | \Xi — Xj\ \ is the distance between the ith and j th trajectories in 
m-dimensional space and is chosen to be smaller than a prescribed small distance 
r*. — is the Euclidian distance between Xi and Xj after k iterations 

of the dynamics F. For our analyses, we have chosen r* to be equal to 4% of the 
maximum value of | |Xj — Xj| |. In all our calculations, k has been determined from 
the value of t at which the autocorrelation function CL{t) of the time series in a or 
Z falls to 1/e of its peak value. The plots are found to get increasingly compact at 
L=5 as m increases and there is no significant change as we go from m=5 to m=6. 
Hence we conclude that the optimal embedding dimension nio for our time series 
in a{t) at 7 = 150s~^ is 5. 

After obtaining a value of m from the local exponential divergence plots, we have 
plotted 7 versus L (right panel. Fig. 4) for m = 3 to 6. As m increases, 7 is found 
to decrease, and there is no significant change in the 7 versus L plot on changing 
m from m=mo (=5) to rrio+i (=6). We see that 7 has a minimum at L = 5 for m 
= 5 and hence we have chosen the optimal delay time to be 5 and the optimal 
embedding dimension nio to be 5 for 7=150s~^. 

We have calculated the correlation dimensions v of the attractors in state space 
on which the trajectories asymptotically lie by computing the correlation integral 
in m-dimensional space C{R) = ^'!^^^ E5=i ^(^ - ll^i - ^ill) [20] for m=2 
to 6. In the expression for C(R), H(x) is the Heaviside function and — Xj\\ 
is the distance between the pair of points in the m-dimensional embedding 
space. The sum in the above expression gives the number of pairs of trajectories 
separated by a distance less than R. For small R's, C(R) is known to scale as 
C{R) ~ R'^, where the correlation dimension 1/ gives useful information about the 
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FIGURE 4. The left panel shows the loeal exponential divergence plots for typical values of 
m and L, marked as (L,ni) for 7=150s~^. The abscissa is ln(||Xj — Xj\\) and the ordinate is 
ln(||Xi_|_fe — — ^^11). The right panel shows the plots of 7 versus L at m=3 to 6, also 

for 7=150s~^. 7 has a minimum at m=5, L=5. 




'log(R) log(R) 
FIGURE 5. Calculation of the correlation dimensions u from the time series corresponding to 
(a) shear stress u and (b) normal stress Z at 25°C and 7=150s~^, for m=2 to 6 (denoted by 
(m,z/)). The values of v calculated from the slopes of the plots of log(C(R)) versus log(R) are 2.9 
for (7 and 3.0 for Z. 
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FIGURE 6. (a) The correlation dimension u and (b) the Lyapunov exponent A as a function of 
the shear rate. 
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FIGURE 7. Plot of 7 versus k for the time series of the shear stresses a at (a) 7=60 (m=5,L=ll), 
(b)90 (m=5,L=9) and (c)150s"^(m=5,L=5). (d) shows the same for the time series of the normal 
stress Z at 7=150s~^ (m=5,L=8). The Lyapunov exponents which can be estimated from the 
slope of the curves are found to be positive for all shear rates shown. 

local structure of the attractor [20]. The plot of C(R) versus R has a constant 
slope over a given range of R, which does not change on increasing m from mo 
to mo+i and is denoted by v. If z/ < m, then the signal exhibits deterministic 
chaos. Fig. 5 shows the calculation of the correlation dimensions for the time 
series in viscoelastic stress a and normal stress Z at 7=150s^^. The slope of the 
C(R) versus R plot in the plateau region is found to saturate at m=5 in both cases 
yielding correlation dimensions v of 2.9 and 3, respectively (Figs. 5(a)-(b)). The 
satisfaction of the condition v <m, and the value of v >2 points to the existence 
of deterministic chaotic dynamics [20] in the stress relaxation of CTAT. We have 
also calculated the correlation dimensions of the attractors corresponding to the 
dynamics of stress relaxation at 60s~^ and 90s~^. The calculated values of v in the 
two cases are 2.4 and 2.7 respectively, indicating the presence of chaotic dynamics 
at these shear rates. Fig 6(a) shows the calculated values of correlation dimensions 

as a function of the shear rate, v is found to increase monotonically with the 
control parameter (shear rate in our case), similar to that observed in the weakly 
turbulent Couette- Taylor flow exhibited by orange oil [21], where the Rayleigh 
number was the control parameter. 

A chaotic time series is characterised by a positive Lyapunov exponent, which 
describes the divergence of neighbouring trajectories in state space [13]. We have 
calculated the largest Lyapunov exponents for the time series of the viscoelastic 
stresses at 60s~^, 90s~^ and 150s~^ and the normal stress at 150s~^ (Fig 7) using 
the method proposed by Gao and Zheng [22]. This method requires the calculation 



of 7(A;) as a function of k, with \\Xi — < and and Xj chosen such that 
\j — i\ > where cu — {m — 1)L, to exclude tangential motion. A linear fit of 
7(A;) versus k yields a positive slope S and nearly zero intercept, from which the 
Lyapunov exponent A has been extracted by using A = S/Atln2. A is found to 
increase monotonically with shear rate as shown in Figure 6(b). 



IV CONCLUSIONS 



The existence of a finite correlation dimension u and a positive Lyapunov ex- 
ponent A indicates the existence of deterministic chaos in the dynamics of stress 
relaxation in CTAT . This occurs only when the shear rates are high enough and lie 
in the plateau region of the flow curve. Since the volume fraction of CTAT in our 
experiment is very small, we rule out the possibility of an isotropic-nematic phase 
transition as the cause of the observed instabilities. We conclude that the chaotic 
dynamics is a natural consequence of the mechanical instability. The calculations 
described above predict an optimal embedding dimension ruo of 5 for the descrip- 
tion of the dynamics of the system. We have previously estabhshed the presence of 
chaotic dynamics in the stress relaxation in CTAT 1.35wt.%, when the sample was 
subjected to step shear rates lying in the plateau region of the flow curve [23]. How- 
ever, at concentrations higher than 2.5wt.%, it is difficult to obtain a sufficiently 
long time series, as the large normal stresses that develop in the sample [2,8] cause 
its expulsion from the rheometer. 

In order to theoretically model the dynamical aspects of the observed phe- 
nomenon in viscoelastic gels, we need to set up space and time-dependent, coupled, 
nonlinear differential equations in at least 5 dynamical variables. One possible 
choice could be the Johnson-Segalman (J-S) equation [24] which is the solution 
of the Navier-Stokes equation coupled to the continuity equation for isothermal, 
incompressible and viscoelastic fluids. The viscoelastic nature of the polymer is 
accounted for by writing the total stress as the sum of a Newtonian part and a de- 
formation history-dependent viscoelastic part. This model yields a non monotonic 
flow curve and oscillations in the start-up stress in polymeric systems [25]. This 
model, however, has not been analysed to predict the chaotic time-dependence of 
the stresses that is seen in our experiments. We believe that the flow-concentration 
coupling and the dynamics of the mechanical interfaces between shear bands need 
to be incorporated in the J-S model. Further, the reversible breakdown and recom- 
bination of micelles also need to be considered. A model constructed by incorpo- 
rating these additional features, which takes into account the nonlinear coupling 
between the relevant dynamical variables like shear and normal stress, shear rate 
and concentration profiles is likely to exhibit the chaotic behaviour that we observe. 
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